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ABSTRACT 

We discuss a novel sampling theorem on the sphere developed by McEwen & Wiaux recently through an associ- 
ation between the sphere and the torus. To represent a band-limited signal exactly, this new sampling theorem 
requires less than half the number of samples of other equiangular sampling theorems on the sphere, such as 
the canonical Driscoll & Healy sampling theorem. A reduction in the number of samples required to represent 
a band-limited signal on the sphere has important implications for compressive sensing, both in terms of the 
dimensionality and sparsity of signals. We illustrate the impact of this property with an inpainting problem on 
the sphere, where we show superior reconstruction performance when adopting the new sampling theorem. 
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1. INTRODUCTION 

The fast Fourier transform^ (FFT) is arguably the most important and widely used numerical algorithm of our 
era, rendering the frequency content of signals accessible. Moreover, Shannon's sampling theorerrP states that 
all of the information content of a band-limited continuous signal can be captured through a hnite number of 
samples. Typically, standard Fourier analyses are performed in Euclidean space where Shannon's theory holds 
and where FFTs are directly applicable. However, in many applications data are observed on non-Euclidean 
manifolds, such as the sphere. Fourier analysis is performed on the sphere in the basis of spherical harmonics, 
which are the eigenfunctions of the spherical Laplacian operator and form the canonical orthonormal basis on 
the sphere. Sampling theorems and fast algorithms to perform spherical harmonic analyses exist but the held is 
much less mature that its elder Euclidean sibling. 

A novel sampling theorem on the sphere has been developed recently by two of the authors of this article^ 
(hereafter referred to as the MW sampling theorem). From an information theoretic viewpoint, the fundamental 
property of any sampling theorem is the number of samples required to capture all of the information content of 
a band-limited signal. To represent exactly a signal on the sphere band- limited at L, all sampling theorems on 
the sphere require 0(L 2 ) samples. However, the MW sampling theorem requires ~ 2L 2 samples only, less than 
half of the ~ 4L 2 samples of other equiangular sampling theorems on the sphere, such as the canonical Driscoll 
& Healy sampling theorerrP (hereafter referred to as the DH sampling theorem) . Not only is the MW sampling 
theorem of theoretical interest, particularly regarding the information content of signals on the sphere, but it 
also has important practical implications in the emerging field of compressive sensing. 

The theory of compressive sensing states that it is possible to acquire sparse or compressible signals with 
fewer samples than standard sampling theorems would suggest P^l T n these settings, the ratio of the required 
number of measurements to the dimensionality of the signal scales linearly with its sparsityP The more efficient 
sampling of the MW sampling theorem reduces the dimensionality of the signal in the spatial domain, thereby 
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improving the performance of compressive sensing reconstruction on the sphere when compared to alternative 
sampling theorems.^ Furthermore, for sparsity priors defined in the spatial domain, such as signals sparse in 
the magnitude of their gradient, sparsity is also directly related to the sampling of the signal. For this class 
of signals, an additional enhancement in compressive sensing reconstruction performance is thus achieved when 
adopting the MW sampling theorem.^ 

In this article we first review sampling theorems on the sphere in Sec. [5J focussing on the MW and DH 
sampling theorems. In Sec. [3] we discuss the superior performance achieved when solving compressive sensing 
problems on the sphere using the MW sampling theorem, as opposed to the DH sampling theorem. We illustrate 
our arguments with an inpainting problem on the sphere, where we adopt the prior assumption that the signal 
to be recovered is sparse in the magnitude of its gradient. Simulations are performed, verifying our theoretical 
arguments. Finally, concluding remarks are made in Sec. 2] 



2. SAMPLING THEOREMS ON THE SPHERE 

Sampling theorems on the sphere state that all of the information contained in a band-limited signal may be 
represented by a finite set of samples in the spatial domain. On the sphere, unlike Euclidean space, the number 
of samples required in the harmonic and spatial domains differ, with different sampling theorems on the sphere 
requiring a different number of samples in the spatial domain. For an equiangular sampling of the sphere, the 
DH sampling theorem has become the standard, while the MW sampling theorem has emerged only recentlyQ 
The MW sampling theorem achieves a more efficient sampling of the sphere, with a reduction by a factor of two 
in the number of samples required to represent a band-limited signal|j In this section we outline the harmonic 
structure of the sphere in the continuous setting, before reviewing concisely the DH and MW sampling theorems. 

2.1 Harmonic Analysis on the Sphere 

We consider the space of square integrable functions on the sphere L 2 (S 2 ) , with the inner product of /, g £ L 2 (S 2 ) 
defined by 

(f,g)= [ Ml(0,<P)W,<P)9*(e,<P), (1) 
Js 2 

where dO(6*, <p) = sinOdOdip is the usual invariant measure on the sphere and (9, ip) define spherical coordinates 
with colatitude 9 £ [0,7r] and longitude <p £ [0, 2ir). Complex conjugation is denoted by the superscript *. 

The spherical harmonic functions form the canonical orthogonal basis for the space of L 2 (S 2 ) functions on 
the sphere and are defined by 



for natural I £ N and integer m £ Z, \m\ < £, where Pg n {x) are the associated Legendre functions. We 
adopt the Condon- Shortley phase convention, with the (— l) m phase factor included in the definition of the 
associated Legendre functions. The orthogonality and completeness relations for the spherical harmonics read 



t=a m=- 



Y lm {9,v)Yt m (p,<ff) = 6(a»9-aH#)8(<p-<ff) (3) 



*Fast algorithms have been developed to compute the forward and inverse transforms rapidly for both the and 
sampling theorems; these algorithms are essential to facilitate the application of these sampling theorems at high 
band-limits. 

^Gauss-Legendre (GL) quadrature can also be used to construct an efficient sampling theorem on the sphere, with 
~ 2L 2 samples.^ The MW sampling theorem nevertheless requires fewer samples and so remains more efficient, especially 
at low band-limits. Furthermore, it is not so straightforward to define norms describing spatial priors on the GL grid 
since it is not equiangular. Finally, algorithms implementing the GL sampling theorem have been shown to be limited to 
lower band-limits and less accurate than the algorithms implementing the MW sampling theorem.^ Consequently, we do 
not focus on the GL sampling theorem any further in this article. 
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respectively, where Sij is the Kronecker delta symbol and S(x) is the Dirac delta function. 

Due to the orthogonality and completeness of the spherical harmonics, any square integrable function on the 
sphere / G L 2 (S 2 ) may be represented by its spherical harmonic expansion 

oo i 

/(0,¥O=E E ftmY im (9,<p) , (4) 

t=0 m=—t 

where the spherical harmonic coefficients are given by the usual projection onto each basis function: 

ftm = (/, Y tm ) = ( dfi(0, ip) f(9, <p) Y; m (0, <p) . (5) 
Js 2 

Throughout, we consider signals on the sphere band-limited at L, that is signals such that fg m = 0, V£ > L. 
In this case the summation over £ in Eqn. (HJ may be truncated to L — 1. We also adopt the convention that 
harmonic coefficients are defined to be zero for \m\ > t (which enforces the contraint \m\ < I when summations 
are interchanged). 

The forward and inverse spherical harmonic transforms, given by Eqn. ([5]) and Eqn. ((4]) respectively, are exact 
in the continuous setting. A sampling theorem on the sphere states how to sample a band-limited function f(9, ip) 
at a finite number of locations, such that all of the information content of the continuous function is captured. 
Since the frequency domain of a function on the sphere is discrete, the spherical harmonic coefficients describe the 
continuous function exactly. A sampling theorem thus describes how to exactly recover the spherical harmonic 
coefficients jg m of the continuous function from its samples. Consequently, sampling theorems effectively encode 
(often implicitly) an exact quadrature rule for evaluating the integral of a band-limited function on the sphere. 

2.2 Driscoll & Healy Sampling Theorem 

The DH sampling theorem^ gives an explicit quadrature rule for the evaluation of spherical harmonic coefficients: 

2L-1 2Z-1 

ftm -EE « DH ^) Vf) Y L(QU Vv) , (6) 

t=0 p=0 

where the equiangular sample positions are defined by 9t — nt/2L, for t = 0, . . . , 2L — 1, and (p p = np/ ' L, for 
p = 0, . . . , 2L— 1, giving TVdh = (2£ — l)2L + 1 ~ 4i 2 samples on the sphere^ The quadrature used to evaluate 
Eqn. ([5]) is exact for a function band-limited at L, with quadrature weights given implicitly by the solution to 

2i-i 2 

]T cornet) P t (cos6 t ) = — 5 e0l V£<2L. (7) 
t=o 

The quadrature weights satisfying Eqn. (J7J are given by 

m , 2tt ^ sin((2fc + l)9 t ) 

fc=0 

The exactness of the quadrature rule is proved by considering the sampling distribution of Dirac delta 
functions defined by 

2L-1 2L-1 

s(9,p)=^2 ^2 «DH(<?t) <5(cos 9 - cos 9 t ) 5{ip - (p p ) . (9) 
t=o p=0 



*The original DH sampling theorem has been revisitedP yielding an alternative formulation with only very minor 
differences and that also requires ~ 4L 2 samples. 
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For quadrature weights satisfying Eqn. 0, it can be shown that sqq = Vw and S£ m = for < £ < 2L, Vm. 
Consequently, the sampling distribution may be written 

oo £ 

s(9, <p) = l+J2 E Sem Yim & v) ■ ( 10 ) 

e=2L m=-l 

The harmonic coefficients of the product of the original band-limited function and the sampling distribution 
f s = f-s are given by 

2L-1 2L-1 

fim -EE 9dh(^) m, p p ) Y; m (9 t , <p v ) , (11) 

t=0 p=0 

which follows from Eqn. (O. Notice that these harmonic coefficients are given by the quadrature rule specified 
in Eqn. ([6]) and it simply remains to prove that the harmonic coefficients of / s agree with those of / for the 
harmonic range of interest (i.e. for < I < L). Noting Eqn. ([TO]) , we may write f s (9,<p) = f(6,<p) + a(6,(p), 
where 

oo I L-l I' 

a(9,cp)=J2 E «tmY tm (6,<p) E E ( 12 ) 

t=2L m=-l £'=0 m'=-l' 

Since the product of two spherical harmonic functions Yg m (6, ip) Y£> m >(8, ip) can be written as a sum of spherical 
harmonics with minimum degree \£ — f'lJ^the aliasing error a(9 7 ip) contains non-zero harmonic content for t > L 
only. Aliasing is therefore outside of the harmonic range of interest and // TO = fi m for < i < L, \m\ < £, thus 
proving the exact quadrature rule given by Eqn. ©. 

2.3 McEwen & Wiaux Sampling Theorem 

The MW sampling theorerrP follows by a factoring of rotationiP and a periodic extension in colatitude 9, so that 
the orthogonality of the complex exponentials over [0, 2tt) may be exploited. This approach encodes an implicit 
quadrature rule on the sphere, which can then be made explicit. 

The spherical harmonics are related to the Wigner functions through^ 



/ 0/ _r_ 1 

Y tm (9, V ) = ^-^-1^^,9,0), (13) 

where the Wigner functions form the canonical orthogonal basis on the rotation group SO (3). The Wigner 
functions may be decomposed aJ^ 

Di n (a, 13, 7 ) = e~ ima d e mn ((3) e"^ , (14) 

where the rotation group is parameterised by the Euler angles (a, /?, 7), with a € [0, 2tt), (3 £ [0, n] and 7 £ [0, 2n). 
The Fourier series decomposition of the (i-functions is given bjff^l 

i 

<4„(/3)=i"- m E A^A^e 1 "^, (15) 

m' = -l 

with A^ ra = d^ nn (n/2). This expression follows from a factoring of rotationsP The Fourier series representation 
of d e mn ((3) given by Eqn. (|15[) allows one to write the spherical harmonic expansion of f(9, ip) in terms of a Fourier 
series expansion of the function extended appropriately to the two-torus T 2 . Noting Eqn. (fT3|) - Eqn. (IT5|) . the 
forward spherical harmonic transform may be written 



L-l 



fern —^ m \j 4?r E ^-rn'm ^L'O G mm > , (16) 

m'=-(i-l) 

where 

Gw=/ d9sm9G m {9)e- im ' e (17) 
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and 2 

G m (0)= [ * d<pf(6,<p)e- im v. (18) 
Jo 

Since Eqn. (|18l) is simply a Fourier transform, the discrete and continuous orthogonality of the complex expo- 
nentials may be exploited to evaluate this integral exactly by 

2 L-i 

Gm{ t) = —— £ f(9 t ^ p )e~ im ^, (19) 
P =-(i-i) 

where y> p = 2irp/{2L - 1), for p = 0, . . . , 2L - 2, and 9 t = 7r(2t + l)/(2i - 1), for t = 0, . . . , L - 1, giving 
JVmw = (L — 1)(2L — 1) + 1 ~ 2L 2 samples on the sphere. It remains to develop a quadrature rule to evaluate 
Eqn. (TTT|) . This is achieved by extending G m (6) to the domain 6 [0,27r) through the construction 



G, 



G m (0 t ), {0,1,..., L-l} 

(-l) m G m (9 2L -2-t) , te{L,...,2L-2} ' 



so that G m {9 t ) may be expressed by a Fourier series. The factor (—1)™ is required to ensure that the symmetry 
in the domain [0, 2tt) dictated by the inverse transform is preserved. Substituting the Fourier expansion of 
G m (9 t ) into Eqn. ([T7J| yields 

L-l 

G mm , =2n F mm „ w(m" - m') , (20) 

m" = -(L-l) 

where the weights are given by 

{±i7r/2, to' = ±1 

0, to' odd, m' ^ ±1 , (21) 
2/(1 — to' 2 ), to' even 

with 

L-l 



t=-(L-l) 

Since the spherical harmonic coefficients fi m are recovered exactly, all of the information content of the function 
f(6,(p) is captured in the finite set of samples. 

The derivation above effectively gives an implicit quadrature rule for the exact integration of a band-limited 
function on the sphere. This quadrature rule can be written explicitly as^ 

. L-l 2L-2 

/ ^ dfip, ip) m v ) = E ft«w(«i) m, , (22) 

•^s 2 t=0 p=0 

where the quadrature weights are defined by 

?mw(0*) = [«(0 t ) + (1 - St, L -i) v(9 2L -a-t)] , (23) 

and where v(9t) is the inverse discrete Fourier transform of the reflected weights w(-m'): 

L-l 

2L- 1 



<0 t ) = -^— J2 w(-m')e im '^. (24) 



-(L-l) 
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3. COMPRESSIVE SENSING ON THE SPHERE 



Compressive sensing on the sphere has been studied recently for signals sparse in harmonic space or in a redundant 
set of overcomplete dictionaries P^E^ However, many natural signals are sparse in measures defined in the spatial 
domain, such as in the magnitude of their gradient. A more efficient sampling of a band-limited signal on the 
sphere, as afforded by the MW sampling theorem, improves both the dimensionality and sparsity of the signal in 
the spatial domain, which has been shown to improve the quality of compressive sampling reconstruction.^ We 
review this very recent work, discussing the impact of efficient sampling on the sphere in the context of a total 
variation (TV) inpainting problem, after first defining the discrete TV norm on the sphere. 

3.1 TV Norm on the Sphere 

The continuous TV norm on the sphere is defined by 

||/||tv= / dfi|V/|, 
where the magnitude of the gradient of the signal / is given by 

|V/| = 



\ 2 / \ 2 

df\ 1 f df\ 



In practice, however, one must consider the TV norm of the sampled signal, where the samples of f(0,tp) are 
denoted by the concatenated vector x £ C , where N is the number of samples on the sphere of the chosen 
sampling theorem (hereafter harmonic coefficients fg m are also represented by a concatenated vector, denoted 
x £ C ). A discrete TV norm on the sphere is defined by approximating the continuous norm in the context of 
either the DH or MW sampling theorem. The integral of the continuous TV norm can be approximated using 
the quadrature rule corresponding to the sampling theorem on the sphere adopted: 

No-l Np-l 

II/Htv^ Yl E IWI^t), (25) 

t=0 p=0 

where the number of samples in (9, ip), given by Ng and N v respectively, and the quadrature weights q(&t) depend 
on the choice of sampling theorem. If |V/| were band-limited, then Eqn. (j2"5)l would be exact. Although this 
is not likely to be the case, Eqn. (j2"5)l is nevertheless a reasonable approximation of the continuous TV norm. 
The magnitude of the gradient |V/| can be approximated from the samples x using finite differences, to give 
a discrete TV norm on the sphere that approximates the continuous norm closely: ||:e||tv — ||/||tv- Notice 
that the inclusion of the quadrature weights q(9 t ) regularises the sin 6 term that arises from the definition of the 
gradient on the sphere, eliminating numerical instabilities that would otherwise occur. 

3.2 TV Inpainting on the Sphere 

We illustrate the impact of the number of samples of the DH and MW sampling theorems on compressive sensing 
reconstruction with an inpainting problem, where measurements are made in the spatial domain. A test signal 
sparse in its gradient is constructed from a binary Earth map, smoothed to give a signal band-limited at L = 32 
(see Fig.[[](a))l The real inpainting problem y — Qx + n is considered, where M noisy measurements y £ R M of 
the signal on the sphere x £ R N are made. The measurement operator $ £ R MxAr represents a random masking 
of the signal. The noise n £ R M is assumed to be independent and identically distributed Gaussian noise, with 
zero mean. 



"The original Earth topography data are taken from the Earth Gravitational Model (EGM2008) pub- 
licly released by the U.S. National Geospatial-Intelligence Agency (NGA) EGM Development Team. 
These data were downloaded and extracted using the tools available from Frederik Simons' webpage: 
http : //www . princeton . edu/ geosciences/people/simons/ 
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(a) Ground truth (b) DH reconstruction (c) MW reconstruction 

Figure 1. Earth topographic data reconstructed in the harmonic domain for M/L 2 = 1/2 



The TV inpainting problem is first solved directly on the sphere: 



x* = argmin ||£e||tv such that \\y — <1>£C 1 1 2 < £ 
x 



(26) 



where the bound e is related to a residual noise level estimator. By adopting the MW sampling theorem in place 
of the DH sampling theorem, the dimensionality and sparsity of the signal in the spatial domain is optimised. 
However, no sampling theorem on the sphere reaches the optimal number of samples in the spatial domain 
suggested by the L 2 dimensionality of the signal in the harmonic domain. Consequently, the dimensionality of 
the problem is reduced by recovering the harmonic coefficients x directly: 



where A € C NxL represents the inverse spherical harmonic transform; the signal on the sphere is recovered by 
x* = Ax* . For this problem the dimensionality of the signal directly recovered x is identical for both sampling 
theorems, however sparsity in the spatial domain remains superior (i.e. fewer non-zero values) for the MW 
sampling theorem. 

Reconstruction performance is plotted in Fig. [5] when solving the inpainting problem in the spatial and 
harmonic domains, through Eqn. (1261) and Eqn. (1271) respectively, for both sampling theorems (averaged over 
ten simulations of random measurement operators and independent and identically distributed Gaussian noise). 
Strictly speaking, compressed sensing corresponds to the measurement ratio M/L 2 < 1 when considering the 
harmonic representation of the signal. Nevertheless, experiments are extended to M/L 2 ~ 2, corresponding 
to the equivalent of complete sampling on the MW grid. When solving the inpainting problem in the spatial 
domain we see a large improvement in reconstruction quality for the MW sampling theorem when compared to 
the DH sampling theorem. This is due to the enhancement in both dimensionality and sparsity afforded by the 
MW sampling theorem in this setting. When solving the inpainting problem in the harmonic domain, we see a 
considerable improvement in reconstruction quality for each sampling theorem, since the dimensionality of the 
recovered signal is optimal in harmonic space. For harmonic reconstructions, the MW sampling theorem remains 
superior to the DH sampling theorem due to the enhancement in sparsity (but not dimensionality) that it affords 
in this setting. In all cases, the superior performance of the MW sampling theorem is evident. In Fig. [1] example 
reconstructions are shown, where the superior quality of the MW reconstruction is again clear. 



Although compressive sensing states that sparse or compressible signals may be acquired with fewer samples 
than standard sampling theorems would suggest, the sampling theorem adopted nevertheless has an important 
influence on the performance of compressive sensing reconstruction. In Euclidean space, Shannon's sampling 
theorem provides an optimal sampling for regular grids, leading to a unique sampling theorem. On the sphere, 
however, no sampling theorem is optimal, with different sampling theorems requiring a differing number of 
samples. The MW sampling theorerrP has been developed only recently and achieves a more efficient sampling 
of the sphere than alternatives, requiring fewer than half as many samples as the canonical DH sampling theorem^ 
while still capturing all of the information content of a band-limited signal. A reduction by a factor of two in 
the number of samples between the DH and MW sampling theorems has important implications for compressive 



x* = argmin ||Ax||tv such that \\y — $Ax\\ 2 < e , 
x 



(27) 



4. CONCLUSIONS 
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Figure 2. Reconstruction performance for the DH and MW sampling theorems 

sensing on the sphere, both in terms of the dimensionality and sparsity of signals. The more efficient sampling 
of the MW sampling theorem has been shown to enhance the performance of compressed sensing reconstruction 
on the sphere, as illustrated with an inpainting problem.^ 
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